clear all;
close all;
clc;

%User cost computation  
r       =.03 ; %Real rate
delta   =.15 ; %Wilson depreciation rate
pi_k    =0.0 ; %Growth of RD input prices
tau     =.07 ; %Giroud and rauh
pk      =1   ; %Price of intangible capital 

%Fsolve -- can also just solve in closed form...
d_target=.026; % target change in user cost (reduction)
d_user_cost=@(c,tau,c_base,d_target,r)  (r+delta-pi_k)*pk*(1-tau-c_base)./(1-tau)-(r+delta-pi_k)*pk*(1-tau-c)./(1-tau)-d_target;

wd=7;
ln=5;
dayy=1;
font=12; 

tau_grid=[.01:.01:.4];
c_grid  =[0.1:.01:.3];
r_grid  =[.01:.005:.1];
c_base  =.25; %RD tax credit base

for i=1:length(r_grid);
    for j=1:length(tau_grid);
        c_out(i,j) = fsolve(@(c) d_user_cost(c,tau_grid(j),c_base,d_target,r_grid(i)),[.35])
        credit_diff_r(i,j)=c_out(i,j)-c_base;
    end
end


f01=figure('name','Figure 12')
mesh(tau_grid,r_grid,credit_diff_r)
xlabel('Corporate tax rate')
ylabel('Interest rate')
zlabel('\Delta RD tax credit rate')
h=title('Change in RD tax credit rate equivalent to 2.6% lower user cost')
view([-220,30])
myfiguresize=[0,0,wd,ln];
set(h,'fontsize',font,'FontWeight','Normal','fontname','times');
set(gca, 'fontsize',font,'FontWeight','Normal','fontname','times');
set(gcf, 'PaperPosition', myfiguresize);
eval([ 'print(f01,''-r300'',''-dpng'',''Figure_12_' num2str(dayy) '.png'')' ])


